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Abstract 

Jamming in hard-particle packings has been the subject of considerable interest in recent years. 
In a paper by Torquato and Stillinger [J. Phys. Chem. B, 105 (2001)], a classification scheme of 
jammed packings into hierarchical categories of locally, collectively and strictly jammed configura- 
tions has been proposed. They suggest that these jamming categories can be tested using numerical 
algorithms that analyze an equivalent contact network of the packing under applied displacements, 
but leave the design of such algorithms as a future task. In this work we present a rigorous and 
efficient algorithm to assess whether a hard-sphere packing is jammed according to the afformen- 
tioned categories. The algorithm is based on linear programming and is applicable to regular as 
well as random packings of finite size with hard-wall and periodic boundary conditions. If the 
packing is not jammed, the algorithm yields representative multi-particle unjamming motions. We 
have implemented this algorithm and applied it to ordered lattices as well as random packings of 
disks and spheres under periodic boundary conditions. Some representative results for ordered and 
disordered packings are given, but more applications are anticipated for the future. One important 
and interesting result is that the random packings that we tested were strictly jammed in three 
dimensions, but not in two dimensions. Numerous interactive visualization models are provided 
on the authors' webpage. 



Electronic address: torquato@electron.princeton.edu 



I. INTRODUCTION 



Packings of hard particles interacting only with infinite repulsive pairwise forces on con- 
tact are applicable as models of complex manybody systems because repulsive interactions 
are the primary factor in determining their structure. Hard-particle packings are therefore 
widely used as simple models for granular media [1], 0, glasses ||, liquids and other 
random media |J, to mention a few examples. Furthermore, hard-particle packings, and 
especially hard-sphere packings, have inspired mathematicians and been the source of nu- 
merous challenging (many still open) theoretical problems |J . Of particular theoretical and 
practical interest are jammed configurations of hard particles. The statistical physics of 
large jammed systems of hard particles is a very active field of theoretical research. For 
packings of smooth hard spheres a rigorous mathematical foundation can be given to the 
concept of jamming, though such rigor is often lacking in the current physical literature. 

There are still many important and challenging questions open even for the simplest 
type of hard-particle packings, i.e., monodisperse packing of smooth perfectly impenetrable 
spheres. One category of open problems pertains to the enumeration and classification 
of disordered disk and sphere packings, such as the precise identification and quantitative 
description of the maximally random jammed (MRJ) state J7j, which has supplanted the 
ill-defined "random close packed" state. Others pertain to the study of ordered systems 
and finding packing structures with extremal properties, such as the lowest or highest (for 
polydisperse packings) density jammed disk or sphere packings, for the various jamming 
categories described below |9|. Numerical algorithms have long been the primary tool for 
studying random packings quantitatively. In this work we take an important step toward 
future studies aimed at answering the challenging questions posed above by designing tools 
for algorithmic assesment of the jamming category of finite packings. 

In the first part of this paper, we present the conceptual theoretical framework underlying 
this work. Specifically, we review and expand on the hierarchical classification scheme for 
jammed packings into locally, collectively and strictly jammed packings proposed in Ref. 
flO|| . In the second part, we present a randomized linear programming algorithm for finding 
unjamming motions within the approximation of small displacements, focusing on periodic 
boundary conditions. Finally, in the third part we apply the algorithm to monodisperse 
packings under periodic boundary conditions, and present some representative but non- 
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exhaustive results for several periodic ordered lattice packings as well as random packings 
obtained via the Lubachevsky-Stillinger packing algorithm JTT . 



Ideas used here are drawn heavily from literature in the mathematical community (fll^, [13 



15 ], etc.), and these have only recently percolated into the granular materials community 



([[H)|, |TE| , [TT| ], etc.). A separate paper will attempt to give a unified and more rigorous 
presentation of the mathematical ideas underlying the concepts of stability, rigidity and 
jamming in sphere packings. Nonetheless, some mathematical preliminaries are given here, 
a considerable portion of which is in the form of footnotes. 

A. Jamming in Hard-Sphere Packings 

The term jamming has been used often with ambiguity, even though the physical intuition 
behind it is strong: It imparts a feeling of being frozen in a given configuration. Two main 
approaches can be taken to defining jamming, kinematic and static. In the kinematic 
approach, one considers the motion of particles away from their current positions, and this 
approach is for example relevant to the study of flow in granular media 1 . The term jammed 
seems most appropriate here. In the static approach, one considers the mechanical properties 
of the packing and its ability to resist external forces 2 . The term rigid is often used among 
physicists in relation to such considerations. 

However, due to the correspondence between kinematic and static properties, i.e. strains 
and stresses, these two different views are largely equivalent. A more thorough discussion of 



this duality is delayed to a later work |18|], but is touched upon in section |VB[ In this paper 
we largely adopt a kinematic approach, but the reader should bear in mind the inherent ties 
to static approaches. 

B. Three Jamming Categories 



First we repeat, with slight modifications as in Ref. the definitions of several hier- 



archical jamming categories as taken from Ref. []10| |, and then make them mathematically 



specific and rigorous for several different types of sphere packings: 



1 In particular, the cessation of flow as jamming is approached. 

2 In particular, the infinite elastic moduli near jamming. 



3 



A finite system of spheres is: 

Locally jammed Each particle in the system is locally trapped by its neighbors, i.e., it 
cannot be translated while fixing the positions of all other particles. 

Collectively jammed Any locally jammed configuration in which no subset of particles 
can simultaneously be displaced so that its members move out of contact with one 
another and with the remainder set. An equivalent definition is to ask that all finite 
subsets of particles be trapped by their neighbors. 

Strictly jammed Any collectively jammed configuration that disallows all globally uniform 
volume-nonincreasing deformations of the system boundary. Note the similarity with 
collective jamming but with the additional provision of a deforming boundary. This 
difference and the physical motivations behind it should become clearer in section 

|vg. 

Observe that these are ordered hierarchically, with local being a prerequisite for collective 
and similarly collective being a prerequisite for strict jamming. 

The precise meaning of some of the concepts used in these definitions, such as unjamming, 
boundary and its deformation, depends on the type of packing one considers and on the par- 
ticular problem at hand, and we next make these specializations more rigorous for specific 
types of packings of interest. Moreover, we should point out that these do not exhaust all 
possibilities and various intricacies can arise, especially when considering infinite packings, 
to be discussed further in Ref. fT8f . It should be mentioned in passing that jammed ran- 
dom particle packings produced experimentally or in simulations typically contain a small 
population of "rattlers", i.e., particles trapped in a cage of jammed neighbours but free to 
move within the cage. For present purposes we shall assume that these have been removed 
before considering the (possibly) jammed remainder. 

C. Unjamming Motions 

Before discussing the algorithms and related issues, it is important to specify exactly what 
we mean by unjamming. First, it is helpful to define some terms. A sphere packing is a 
collection of spheres in Euclidian rf-dimensional space such that the interiors of no two 
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spheres overlap. Here we focus on monodisperse systems, where all spheres have diameter 
D = 2R, but most of the results generalize to polydisperse systems as well. A packing -P(R) 
of N spheres is characterized by the arrangement of the sphere centers 3 R = (r l5 . . . , r/v), 
called configuration: 

P(R) = {r; e U d , i = 1, . . . ,N : \\n - > D Vj ^ i} 

The natural physical definition 4 of what it means to unjam a sphere packing is provided 
by looking at ways to move the spheres from their current configuration. This leads us to the 
definition: An unjamming motion AR(t), where t is a time-like parameter, t e [0, 1], is a 
continuous displacement of the spheres from their current position along the path R+AR(£), 
starting from the current configuration, AR(0) = 0, and ending at the final configuration 
R+ AR(1), while observing all relevant constraints along the way 5 , such that some of the 
contacting spheres lose contact with each other for t > 0. If such an unjamming motion 
does not exist, we say that the packing is jammed 6 . It can be shown (see references in Ref. 



T4| ) that an equivalent definition 7 is to say that a packing is jammed if it is isolated in the 



allowed configuration space, i.e., there is no valid packing within some (possibly small) finite 
region around R that are not equivalent 8 to -P(R). We mention this because it is important 
to understand that although we use a kinematic description based on motion through time, 
which imparts a feeling of dynamics to most physicists, a perfectly equivalent static view is 
possible. Therefore, henceforth special consideration will be given to the final displacement 9 
AR(1), so that we will most often just write AR = AR(1). 

It is also useful to know that we need to only consider analytic AR(t), which gives us the 
ability to focus only on derivatives of AR(i). Furthermore, it is a simple yet fundamental 
fact 10 that we only need to consider first derivatives 11 V = |AR(t), which can be thought 



3 Capitalized bold letters will be used to denote dTV-dimensional vectors which correspond to the d- 
dimensional vectors of all N of the particles. 

4 This is the second definition (definition b) in section 2.1 of Ref. |]l4j| . 

5 This means that impenetrability and any other particular (boundary) conditions must be observed, i.e. 
P(R + AR(i)) is a valid packing for all t e [0, 1]. 

6 Of course by changing the (boundary) constraints we get different categories of jamming, such as local, 
collective and strict. 

7 This the first definition (definition a) in section 2.1 of Ref. fl4|| . 

8 Two packings are equivalent if there is a distance-preserving (rigid-body) motion (map), such as a uniform 
translation or a rigid rotation, that takes one packing to the other. 

9 This will be important when discussing random packings with finite (but small) interparticle gaps. 

At the second derivative spheres cannot interpenetrate, and in fact strictly move away from each other. 
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of as "velocities", and then simply move the spheres in the directions V = (vi, . . . , V2) to 
obtain an unjamming motion AR(i) = Vi. That is, a sphere packing that is not jammed 
can be unjammed by giving the spheres velocities V such that no two contacting spheres i 
and j, || rj — rj|| = D, have a relative speed Vij toward each other 12 : 

%j = (Vi - Vj) T < (1) 

where Ujj = j^ J ~^ is the unit vector connecting the two spheres 13 . Of course, some special 
and trivial cases like rigid body translations (V = constant) need to be excluded since they 
do not really change the configuration of the system. 

In this paper we will plot unjamming motions as "velocity" fields, and occasionally supple- 
ment such illustrations with a sequence of frames from t — to t — 1 showing the unjamming 
process. Note that the lengths of the vectors in the velocity fields have been scaled to aid in 
better visualization 14 . For the sake of clear visualization, only two-dimensional examples will 
be used. But all of the techniques described here are fully applicable to three-dimensional 
packings as well. Interactive Virtual Reality Modeling Language (VRML) animations which 
are very useful in getting an intuitive feeling for unjamming mechanisms in sphere packings 
can be viewed from Windows platforms on our webpage (see Ref. [21]). 



II. BOUNDARY CONDITIONS 



As previously mentioned, the boundary conditions imposed on a given packing are very 
important, especially in the case of strict jamming. Here we consider two main types of 
packings, depending on the boundary conditions. 

Hard-wall boundaries The packing -P(R) is placed in an impenetrable concave 15 hard- 
wall container K. Figure H shows that the honeycomb lattice can be unjammed inside 



The formal proof and a more detailed discussion will be given in Ref. 18 . 

11 The formal statement is that a packing is rigid if and only if it is infinitesimally rigid, see Refs. jl2[ [uj| . 

12 This is the third definition (definition c) in section 2.1 of Ref. Jl4| . 

13 The sign notation may be a bit unorthodox but is taken from Ref. |2(| . 

14 Since this paper deals only with small displacements in an approximate linearized fashion, one should do 
a line-search in t along R + AR(i) to check when the first collision occurs. The VRML animations shown 
at the webpage of Ref. [21] also depict arbitrarily scaled displacements for visualization purposes so that 
collisions might happen before t = 1. 

15 A container is concave if its boundary is concave at every circular point. All convex polyhedral sets make 
a concave container, but a spherical container is not concave. The details of the definition and the reasons 
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a certain hard-wall container. Often we will make an effective container out of Nf 
fixed spheres whose positions cannot change. This is because it is often hard to fit a 
packing into a simple container such as a square box, while it is easy to surround it with 
other fixed spheres, particularly if a periodic lattice is used to generate the packing. 
Specifically, one can take a finite sub-packing of an infinite periodic packing and freeze 
the rest of the spheres, thus effectively making a container for the sub-packing. An 
example is depicted in Fig. ||. 

Periodic boundaries Periodic boundary conditions are often used to emulate infinite sys- 
tems, and they fit the algorithmic framework of this work very nicely. A periodic 
packing P(R) is generated by a replicating a finite generating packing -P(R) on a lat- 
tice A = {Ai, . . . , Xd}, where Aj are linearly independent lattice vectors and d is the 
spatial dimensionality. So the positions of the spheres are generated by, 

r l{n c ) = + ^ n c an d n c is integer, n c G Z d 

where we think of A as a matrix having the lattice vectors as columns 16 and n c is 
the number of replications of the unit cell along each basis direction 17 . The sphere 
i (n c ) is the familiar image sphere of the original sphere i = i (0), and of course for the 
impenetrability condition only the nearest image matters. For true periodic boundary 
conditions, we wrap the infinite periodic packing -P(R) around a flat torus 18 , i.e. we 
ask that whatever happens to a sphere i also happens to all of the image spheres i (n c ) , 
with the additional provision that the lattice may also change by AA: 

Ar ?(n c) = A *i + (AA)n c (2) 
why the container needs to be concave to make some of the theorems used here work are given in Ref. 

s 

16 The matrix A has d 2 elements. 

17 This can also be viewed as a simple way of generating an infinite packing, and one can analyze the 
resulting infinite packing, called a cover of as discussed in Ref. fl^| , however, infinite packings are 
mathematically delicate and will be discussed in Ref. p8|. 

18 We mean a topological torus, not a geometrical one, i.e., a torus where distances are still measured as in 
flat Eucledian space but which has the topology of a curved torus. 
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A. Using Simple Lattices to Generate Packings 

Simple familiar lattices 19 such as the triangular, honeycomb, Kagome and square in two 
dimensions, or the simple cubic (SC), body-centered cubic (BCC), face-centered cubic (FCC) 
and hexagonal-close packed (HCP) in three dimensions, can be used to create a (possibly 
large) packing taking a subsystem of size N c unit cells along each dimension from the 
infinite lattice packing. The properties of the resulting system can be studied with the tools 
developed here, provided that we restrict ourselves to finite N c . Moreover, it is important to 
specify which lattice vectors are to be used. We will usually take them be primitive vectors, 
but sometimes it will be more convenient to use conventional ones (especially for variations 
on the cubic lattice). 

For hard-wall boundary conditions, we can take an infinite packing generated by these 
simple lattices and then freeze all but the spheres inside the window of N c unit cells, thus 
effectively obtaining a hard- wall container. Figure || illustrates an unjamming motion for 
the honeycomb lattice under these conditions. 

For periodic boundary conditions, the generator P(R) can itself be generated using a 
simple lattice 20 . This is not only a convenient way to generate simple finite periodic packings, 
but it is in general what is meant when one asks, for example, to analyze the jamming 
properties of the Kagome lattice under periodic or hard- wall boundary conditions. Figure |3| 
shows a periodic unjamming motion for the Kagome lattice. Notice though that the jamming 
properties one finds depend on how many neighboring unit cells N c are used as the "base" 
region (i.e., the generating packing), and therefore, we will usually specify this number 
explicitly. Some properties are independent of N c , and tailored mathematical analysis can 
be used to show this [[0|, p2| . For example, if we find a periodic unjamming motion for 
N c = (1, 2, . . . ), then we can be sure that for any N c that is an integer multiple of this 
window we can use the same unjamming motion. Correspondingly, if we show that the 
system with N c = (2 ■ 3 ■ 5, 2 ■ 3, . . . ) is jammed, then so must be any system whose size can 
be factored into the same prime numbers along each dimension. We will not consider these 
issues in detail here, but rather focus on algorithmic approaches tailored for finite and fixed 

19 We mean lattices with a basis, to be more precise. 

20 In this case the lattice A is a sub-lattice of the underlying (primitive) lattice A, i.e., X a — [N c ] is an 
integer multiple of the vector of the underlying lattice A Q , a = 1, . . . , d. 
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systems (i.e., N c is fixed and finite), and postpone the rest of the discussion to Ref. 



18 



III. LINEAR PROGRAMMING ALGORITHM TO TEST FOR JAMMING 



Given a sphere packing, we would often like to test whether it is jammed according to 
each of the categories given above, and if it is not, find one or several unjamming motions 
AR(i). We now describe a simple algorithm to do this that is exact for gap-less packings, i.e., 
packings where neighboring spheres touch exactly, and for which the definitions given earlier 
apply directly. However, in practice, we would also like to be able to study packings with 
small gaps, such as produced by various heuristic compression schemes like the Lubachevsky- 
Stillinger algorithm fLlf . In this case the meaning of unjamming needs to be modified so 
as to fit physical intuition. We do this in such a way as to also maintain the applicability 
of our efficient randomized linear programming algorithm using what Roux ||16[ calls the 
approximation of small displacements (ASD). 



A. Approximation of Small Displacements 

As already explained, an unjamming motion for a sphere packing can be obtained by 
giving the spheres suitable velocities, such that neighboring spheres do not approach each 
other. Here we focus on the case when AR(t) = Vt + 0(t 2 ) are small finite displacements 
from the current configuration. Therefore, we will drop the time designation and just use 
AR for the displacements from the current configuration R to the new configuration R = 
R+ AR. 

In this ASD approximation, we can linearize the impenetrability constraints by expanding 
to first order in AR, 

\\ri-rjW = ||(r,-r J ) + (Ar l -Ar J )|| >D (3) 
to get the condition for the existence of a feasible displacement AR, 

(Ar, - Ar ; j 1 u, ; < Al itj for all {i, j} (4) 
where {i, j} represents a potential contact between nearby spheres i and j, AZ,j = ||rj — r, || — 
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D is the interparticle gap 21 , and Ujj = ||r J Zr'|| is the unit vector along the direction of the 
contact {i, j}. For a gap-less packing, we have Al = and the condition (f|) reduces to 
(H), which as we explained along with the conditions AR 7^ and AR 7^ const, is an exact 
condition for the existence of an unjamming motion AR. For packings with finite but small 
gaps though, condition (|4]) is only a first-order approximation. Notice that we only need 
to consider potential contacts {i,j} between nearby, and not all pairs of spheres 22 . The 
complicated issue of how well the ASD approximation works when the gaps are not small 
enough is illustrated in Fig. f|. 

By putting the U^j S clS columns in a matrix of dimension [Nd x N e ], where N e is the 
number of contacts in the contact network, we get the important rigidity matrix 23 of the 
packing A. This matrix is sparse and has two blocks of d non-zero entries in the column 
corresponding to the particle contact {i,j}, namely, in the block row corresponding to 
particle i and — Uy in the block row corresponding to particle j. Represented schematically: 

{iJ} 
i 

-Uij- 

For example, for the four-disk packing shown in Fig. |], and with the numbering of the disks 

21 Called interstice in Ref. Jig . 

22 That is we only consider a contact if 

Iki-r^l <(l + 5)D 

where S is some tolerance for how large we are willing to allow the representative ||Ar|| to be. The 
larger this tolerance, the more possible particle contacts we will add to our constraints, and thus the 
more computational effort we need. Also, the ASD approximation becomes poorer as we allow larger 
displacements. But choosing a very small tolerance makes it impossible to treat systems with moderately 
large interparticle gaps (say of the order of 6 = 10%). 

23 This is in fact the negative transpose of what is usually taken to be the rigidity matrix, and is chosen to 
fit the notation in Ref. p0| . 




10 



depicted in Fig. |5L we have the following rigidity matrix: 





En 


-^13 




D 1 


Ul2 


U13 




D 2 


-U12 






D 3 




-U13 













Using this matrix, we can rewrite the linearized impenetrability constraints as a simple 
system of linear inequality constraints: 



A T AR < Al 



(5) 



The set of contacts {i,j} that we include in (^) form the contact network of the 
packing, and they correspond to a subclass of the class of fascinating objects called tensegrity 
frameworks, namely strut frameworks (see Ref. for details, and also || for a treatment 
of more general packings). Figure ^ shows a small random packing with relatively large gaps 
and the associated contact network. 



1. Boundary conditions 

Handling different boundary conditions within the above formulation is easy. For ex- 
ample, for usual periodic conditions, handling the boundaries merely amounts to adding a 
few columns to the rigidity matrix A with u.^, \ = n 3(nc) — tt for all images j(n c ) which 

II j("c) Ml 

have contacts with one of the original spheres i. These columns correspond to the periodic 
contacts wrapping the packing around the torus. 

For hard-wall boundaries, we would add a potential particle contact to the contact net- 
work from each sphere close to a wall to the closest point on the wall, and fix this endpoint. 
Such fixed points of contact and fixed spheres 24 j are simply handled by transferring the 
corresponding term Arjujj to the right-hand side of the constraints in ([5]). 

Such nodes are called fixed nodes in tensegrity terminology. 
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B. Finding Unjamming Motions 

We are now ready to explain how one can find unjamming motions for a given packing, 
if such exist. But first, we need to refine our definition of an unjamming motion to allow for 
the study of packings with finite gaps. 

1. Unjamming Motions Revisited: Dealing with Finite Gaps 

The problem with packings with small gaps is that according to our previous definition 
of an unjamming motion in section [K], such packings will never be jammed, since it will 
always be possible to move some of the spheres at distances comparable to the sizes of the 
interparticle gaps so that they lose contact with all or most of their neighbors. Clearly we 
wish to only consider the possibility of displacing some of the spheres such that considerable 
gaps appear, larger then some threshold value A/i arge 3> AZ, where Al is a measure of 
the magnitude of the interparticle gaps. Therefore, we have the modified definition: An 
unjamming motion AR(t), t G [0,1], is a continuous displacement of the spheres from 
their current position along the path R + AR(i), starting from the current configuration, 
AR(0) = 0, and observing all relevant constraints along the way, such that some of the 
spheres lose contact 25 with each other in the final configuration R + AR(1) by more then 
a given AZi arge . Again, we exclude any trivial rigid-body type motions such as a uniform 
translation of the spheres from consideration. 

The problem of whether such an unjamming motion exists and how to find one if it does 
is mathematically very interesting if we take the case AZi arge 3> D. But this problem is 
extremely complex due to high non-linearity of the impenetrability constraints and we will 
make no attempt to solve it. It also is not clear that this would be of interest to physicists. 
Instead, we focus our attention on the case when AZi arge is small enough to apply with a 
reasonable degree of accuracy the approximation of small displacements, but large enough 
compared to the interparticle gaps so that the exact value is really irrelevant. 

Under the ASD, we need only worry about linear displacements from the current config- 
uration, AR(£) = Vi, and so we can focus on AR = AR(1). Thus, finding an unjamming 

25 Alternatively, we could ask that some sphere displace by more than Ari arg0 ^> Al. 
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motion 26 simply reduces to the problem of feasibility of a linear system of inequalities, 



A T AR < Al for impenetrability (6) 

> A/i ar g e 3> Al for unjaniHiing (7) 

where we can exclude trivial displacements such as uniform translations by adding additional 
constraints (e.g., demanding that the centroid remains fixed, ^ Ar; = 0). 



3 {i, j} such that 



ar).. 



2. Randomized Linear Programming (LP) Algorithm 

The question of whether a packing is jammed, i.e., whether the system ([F]) is feasible, can 
be answered rigorously 27 by using standard linear programming techniques 28 . If a packing 
is jammed, then this is enough. But for packings which are not jammed, it is really more 
useful to obtain a representative collection of unjamming motions. A random collection of 
such unjamming motions is most interesting, and can be obtained easily by solving several 
linear programs with a random cost vector. 

We adopt such a randomized LP algorithm to testing unjamming [i.e. studying (0)], 



26 Since here we are really focusing only on AR = AR(1), we should change terminology and call Ar an 
unjamming displacement. However, to emphasize the fact that there is a way to continuously move 
the spheres to achieve this displacement, and also that we can always scale down the magnitude of an 
allowed AR arbitrarily, we continue to say unjamming motion. The term "displacement" will be more 
appropriate in section VB. 

27 The gap-less case Al = 0, AZi arge — > + , can be studied rigorously. When gaps are present, of course, the 



'A T AR) r 



gc 

3> Al is mathematically ambiguous and also the ASD approximation becomes 



condition 
inexact. 

28 Solve the following linear program aimed at maximizing the sum of the (positive) gap dilations 
(Al-A^AR).., 

min AR (A T AR), j = min (Ae) T AR (8) 

such that A T AR < Al (9) 

where e is the unit vector, and then look at the magnitudes of the gap dilations (these may be unbounded 
of course) and decide if they are large enough to consider the solution an unjamming motion. Otherwise 
the packing is jammed. Notice that this will usually produce a single unjamming motion, which we have 
found to be rather uninteresting for lattice packings in the sense that it is extremely dependent upon N c . 
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namely, we solve several instances of the following LP in the displacement formulation: 



niaxAR,b T AR for virtual work (10) 
such that A T AR < Al for impenetrability (11) 
|AR| < A_R max for boundedness (12) 

for random loads 29 b, where Ai? max Al is used to prevent unbounded solutions and thus 
improve numerical behavior 30 . Trivial solutions, such as uniform translations of the packing 
AR = const, for periodic boundary conditions, can be eliminated a posteriori, for example 
by reducing AR to zero mean displacement 31 , or added as extra constraints in (|TT1) . We will 
discuss numerical techniques to solve (|TTJ) shortly. 

We then treat any solution AR to ( |TT|) with components significantly larger then Al as 
an unjamming motion. For each b, if we fail to find an unjamming motion, we apply — b 
as a loading also 32 . In our tests we usually set AR max ~ 10 D and treated any displacement 
where some gap dilations (A T AR) . . ~ D as unjamming motions, but as should be clear 



29 The physical interpretation of b as an external load will be elucidated in section VB . 

30 Even for jammed packings, unless b is chosen carefully, the solution of ( pi] ) will be unbounded due to the 
existence of trivial motions such as uniform translations. This will become clear once duality is discussed, 
but mathematically b needs to be in the null-space of A, which usually means it needs to have zero total 
sum and total torque (see chapter 15 in Ref. f20f). 

31 The choice of b also affects the appearance of these trivial solutions in the optimal solution, which is non 
unique, as explained in footnote [50|. 

32 The linearized impenetrability constraints A T AR < Al define a polyhedral set Par of feasible displace- 
ments. Every such polyhedron consists of a finite piece T^R) ^ e convex nu ^ °f extreme points, and 
possibly an unbounded piece Car, a finitely generated polyhedral cone. In some cases this cone will be 
empty (i.e. Car = {0}), but in others it will not, as can be seen in Fig. |J. A mathematically very well 
defined formulation is to take any ray in the cone Car as an unjamming motion, and exclude others, 
however, as Fig. |] shows, the elongated corners of this polyhedron are in fact very likely to be unbounded 
in the true non-linear feasible set of displacements, so we prefer to take any "long" direction in "Par as 
an unjamming motion. 

We note that the randomized LP algorithm proposed here strictly answers the question of whether the 
polyhedral set of feasible displacements contains an unbounded ray just by applying two (nonzero) loads 
b and — b. This is because an attempt to find such a ray will be unsuccessful only if — h G Carj where 
C AR is the conjugate cone of C AR , and in this case b ^ C AR , so that using the load — b will find a ray 
if such a ray exists. Also, we note that one cannot hope to fully characterize the cone of first-order 
unjamming motions C AR (i.e. find its convex hull of generating rays), as this is known to be an NP 
complete problem related to the full enumeration of the vertices of a polyhedron. Our randomized 
approach essentially finds a few sample rays in C AR . 
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from the discussion in footnote 0, the method is not very sensitive to the exact values. One 

can then save these individual unjamming motions and visualize them, or try to combine 
several of these into one "most interesting" unjamming motion 33 , as we have done to obtain 
the figures in this paper. 

We stress that despite its randomized character, this algorithm is almost rigorous when 
used as a test of jamming, in the sense that it is strictly rigorous for gap- less packings and 
also for packings with small gaps as explained in more detail in footnote 

IV. TESTING FOR LOCAL, COLLECTIVE AND STRICT JAMMING: PERI- 
ODIC BOUNDARY CONDITIONS 

A. Local Jamming 

Recall that the condition for a packing to be locally jammed is that each particle be fixed 
by its neighbors. This is easy to check. Namely, each sphere has to have at least d+1 contacts 
with neighboring spheres, not all in the same <i-dimensional hemisphere. This is easy to test 
in any dimension by solving a small linear program, and in two and three dimensions one 
can use more elementary geometric constructions. 

We prefer the LP approach because it is in the spirit of this work and because of its 
dimensional independence, and so we present it here. Take a given sphere i and its set of 
contacts {uj*}, and put these as rows in a matrix Af. Then solve the local portion of © 
in footnote ^8] (using the simplex algorithm): 



which will have an unbounded solution if the sphere i is not locally jammed, as illustrated 
in Fig. [§ 

Of course we can define higher orders of local jamming by asking that each set of n 
spheres be fixed by its neighbors, called n-stability in Ref. [|l^]. However, for n > 1 it 

33 We believe motions in which as many of the spheres move as possible are most useful since then one can 
see multiple unjamming "mechanisms" with one visualization. Therefore, we make a convex combination 
of several unjamming motions AR(b Q ) obtained from several different random loads b Q to obtain one 
such unjamming motion. 




(13) 
(14) 



such that Af Arj < Al^ 
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becomes combinatorially too difficult to check for this. Computationally, we have found 
testing for local jamming using (O) to be quite efficient and simple. 



B. Collective Jamming 

The randomized LP algorithm was designed to test for collective jamming in large pack- 
ings, and in this case the linear program ([□]) that needs to be solved is very large and sparse. 
Notice that boundary conditions are only involved when making the list of contacts in the 
contact network and deciding if certain spheres or contact points are fixed. In the case of 
periodic boundary conditions, we simply add the usual contacts between original spheres 
near the boundary of the unit cell and any nearby image spheres. 

We have implemented an efficient numerical solution of ( |TTD [and also (0)], using the 
primal-dual interior-point algorithm in the LOQO optimization library (see Ref. P5| ). Illus- 
trations of results obtained using this implementation are given throughout this paper. 

We would like to stress that primal-dual interior-point algorithms are very well suited 
for problems of this type, and should also be intuitive to physicists since in essence they 
solve a sequence of easier problems in which the perfectly rigid inter-sphere contacts are 
replaced by stiff (but still deformable) nonlinear (logarithmic) springs, carefully numerically 
taking the limit of infinitely stiff springs. Physicists have often used similar heuristically 
designed schemes and hand-tuned them, and even suggested that standard optimization 
algorithms are not practical (for example, in Ref. |16|]). We would like to dispel such beliefs 



and stress the importance of using robust and highly efficient software developed by applied 
mathematicians around the world, such as Ref. |p3|| , becoming increasingly more available. 



Not only are the algorithms implemented theoretically well-analyzed, but they are tested 
on a variety of cases and often contain several alternative implementations of computation- 
intensive sections targeting different types of problems. Choice of the correct algorithm and 
the details are often complex, but well worth the effort. 

Nonetheless, for three-dimensional problems the available implementations of interior- 
point algorithms based on direct linear solvers are too memory demanding and inefficient. 
Tuned implementations based on conjugate- gradient iterative solvers are needed. We plan 
to develop efficient parallel algorithms suited for these types of problems and make them 
publicly available in the very near future. 
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C. Strict Jamming 



To extend the notion of collective jamming to strict jamming we introduced deformations 
of the boundary. In the case of periodic packings, the lattice A is the boundary. Therefore, 
the only difference with collective jamming is that we will now allow the lattice to change 
while the spheres move, i.e., A A ^ in (0). The lattice deformations A A will also 
become unknowns in (|TID , but since these too enter linearly in (0), we still get a linear 
program, only with coefficient matrix A augmented with new (denser) rows in the columns 
corresponding to contacts across the periodic boundary. The actual implementation now 
requires more care and bookkeeping, but the conceptual changes should be clear, and the 
randomized LP algorithm remains applicable. 

Obviously, we cannot allow the volume of the unit cell to enlarge, since the unit cell is 
in a sense the container holding the packing together. Therefore, we only consider volume- 
non-increasing continuous lattice deformations AA(t): 



det 



A = A + AA(t) < det A fort >0 (15) 



Through a relatively simple mathematical analysis to be presented in Ref. ||18|| , it can 



be shown that the principles that applied to unjamming motions of the spheres AR still 
remain valid even when we extend the notion of an unjamming motion to include a deforming 
lattice 34 . That is, we can still only focus on linear motions AA(t) = Wt, W = const, and the 
final small deformations A A = AA(1), and need to consider only first-order linearizations 
of the impenetrability (|3|) and non-expansion ( |T5| ) nonlinear constraints 35 . 
The linearized version of (|T5|) is: 

Tr[(AA)A" 1 ] < (16) 

and this is just one extra linear constraint to be added to the linear program (|ll|). An extra 
condition which needs to be added is that (AA)A -1 be symmetric, which is also an added 
linear constraint, 

(AA)A- 1 = e where e = e T (17) 



34 That is, we now think of [AR(i), AA(t)] as an unjamming motion. 

35 The condition (I!?]) needs to hold also. 
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where we add e as an unknown in the randomized LP algorithm. This condition in fact does 
nothing more then eliminate trivial rotations 36 of the lattice 37 . 

The motivation for the category of strict jamming and its above interpretation in the 
periodic case should be clear: Changing the lattice in a volume non-increasing way models 
macroscopic non-tensile strain and is therefore of great relevance to studying the macro- 
scopic mechanical properties of random packings (see Ref. |19|| ). In fact, e = (AA)A -1 can 
be interpreted as the "macroscopic" strain-tensor, which explains why it is symmetric and 
also trace-free for shear deformations. We also again point out that strict jamming is (signif- 
icantly) stronger then collective jamming for periodic boundary conditions, particularly in 
two-dimensional packings 38 . This point is illustrated in Fig. |7], which shows an unjamming 
motion involving a deformation of the lattice, even though this lattice packing is collectively 
jammed. 



D. Shrink- And-Bump Heuristic 



The following heuristic test for collective jamming has been suggested in Ref. [IT|: Shrink 
the particles by a small amount 5 and then start the Lubachevsky-Stillinger molecular dy- 
namics algorithm with random velocities, and see if the system gets unjammed. One would 
also slowly enlarge the particles back to their original size while they bump around, so 
as to allow finite termination of this test (within numerical accuracies). We call this the 
shrink- and-bump heuristic. The idea is that the vector of velocities takes on random values 
in velocity space and if there is a direction of unjamming, it will be found with a high 
probability and the system will unjam 39 . Animations of this process can be found at Ref. 
[21]. 

This kind of heuristic has the advantage of being very simple and thus easy to implement 



36 Rotations of the lattice turn out to increase the unit-cell volume at the second-order derivative, even 
though they are volume-preserving up to first order. 

37 But one should still deal with trivial motions a posteriori with some care in certain pathological cases. 

38 This point will be elaborated in Ref. @. 

39 The theory presented here suggests that if a packing is indeed not collectively jammed and has a relatively 
large cone of unjamming motions (see footnote |3^), it can be unjammed using this type of heuristic with 
high probability. However, notice that this cone can in principle be very small, so finding a ray in it 
may be a low-probability occurrence. For packings with finite gaps, though, the heuristic incorporates 
nonlinear effects, which is an advantage. 
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and use, and it is also very efficient. The real problem is not so much its indeterminacy, but 
its strong dependence on the exact value of 5. For example, animations showing how the 
Kagome lattice inside a container made of fixed spheres (as in Fig. |2|) can be unjammed 
with a large-enough 5, even though it is actually collectively jammed under these boundary 
conditions, can be found at Ref. [21]. In fact, many jammed large packings will appear 
unstable under this kind of test, as motivated with the notion of uniform stability, defined 



in Ref. [13] and elaborated on in Ref. 18 



V. ADDITIONAL APPLICATIONS 



A. Compressing Packings using Linear Programming 

In this work we emphasize the utility of the randomized linear programming algorithm 
as a testing tool for jamming, and also for finding representative unjamming motions. The 
unjamming motions one finds can be used inside compression algorithms. For example, the 
Lubachevsky-Stillinger algorithm may sometimes get stuck in a particular configuration even 
though the configuration is not collectively or strictly jammed (particularly in two dimen- 
sions, as explained later). The unjamming motion obtained from the linear programming 
algorithm can then be used to continue the compression 40 . 

Even more can be done with linear programming in this direction. For example, one can 
ask the question of whether there is an unjamming motion AR in which all sphere contacts 
are lost 41 . This can be done for example by solving the LP, 

maxAR l£ a (18) 
such that A T AR < -aD (19) 
< a < 1 for boundedness (20) 

By displacing the spheres by such a AR we would create a "cushion" of free space around 
each sphere, so that we can actually increase the radius of the spheres 42 and thus increase 
the density, or equivalently, the packing fraction <p. We believe these kinds of approaches to 

40 For example, one can displace the spheres by AR and restart the compression with random velocities or 
use initial velocities along AR to unjam the packing and continue the simulation. 

41 All contacts will be lost if A T AR < —aD, where a > 0. 

42 The common radius could increase by at least aD. 
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be too inefficient, since solving a large-scale linear program is too expensive to be used iter- 
atively. This is reminiscent of the high-quality but rather expensive compression algorithm 
of Zinchenko (see Ref. p4|j). 

Our future work will focus on using mathematical programming algorithms and rigidity 
theory to design high-quality algorithms for design of packings with target properties, using 
systems of stiff nonlinear springs as an intermediary. 



B. Kinematic/Static Duality 

The subject of kinematic/static duality and its physical meaning and implications are 
discussed at length in Ref. and elsewhere in various degrees of relevance and different 



perpectives fT0| , |12| , |14| , |15| , |17| . Here we only comment on it because of its relevance to the 
randomized LP algorithm for testing jamming, and leave further discussion of this important 
subject to Ref. JT8| . 



The dual 43 of the displacement formulation LP ((TTJ) also has a very physical interpretation 
and it gives us the interparticle repulsive 44 forces f as dual variables, and we call it the force 
formulation LP: 

maxf (Al) T f for virtual work (21) 
such that Af = b for equilibrium (22) 
f < for repulsion only (23) 

The physical interpretation of the objective function in both the displacement (|ll|) and 
force formulations ([22|) is that of (virtual) mechanical work done by the external force load 
b applied to the spheres. These two LP's are of great importance in studying the stress- 
strain behavior of granular materials, as explained in Ref. |16j, and since they are equivalent 
to each other, we can call them the ASD stress-strain LP. 

We wish to emphasize that by using primal-dual interior point algorithms we automati- 
cally get both forces and displacements using the same implementation 45 . We have empha- 
sized the displacement formulation ( |i"Tl) simply because we based our discussion of jamming 



43 Excluding the additional practical safeguard constraint AR < AR max , which is added to avoid unbounded 
trivial or unjamming motions. 

44 We choose a negative sign for repulsive forces here in agreement with mathematical literature 

45 For example, both LOQO and PCx (Q) return both primal and dual solutions to the user. 
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on a kinematic perspective, but a parallel static interpretation can easily be given. For 
example, a random b used in the randomized LP algorithm that finds an unbounded un- 
jamming motion physically corresponds to a load that the packing cannot support, i.e. the 
force formulation LP is (dual) infeasible, implying that the displacement formulation LP is 
(primal) unbounded 46 . The meaning of collective jamming within the ASD in the presence of 
small gaps from a static standpoint now becomes clear: A collectively jammed packing can 
resist (support) any force loading by (as) small (as possible) rearrangements of the spheres, 
in which some of the potential contacts are open and others closed, depending on the loading 
and the interparticle gaps. 

In general the stress-strain LP will be highly degenerate and its primal and/or dual 
solution not unique. However, as Roux points out, the existence of small gaps in random 
packings is very important in this context. Namely, if Al is random and nonzero (however 
small), and b is also random, theorems on the generic character of linear programs (see the 
references in Ref. [20]) can be invoked to guarantee that both the primal and dual solutions 
will be non-degenerate. A non-degenerate solution to (p2|) corresponds to an isostatic force- 
carrying contact network, a fact noted and explained in a great many ways by various 
researchers in the field of granular materials |15|, [IB, |T7). We just mention these points here 
in order to stimulate interest among the physical community in the very relevant results to 
be found in the mathematical programming literature. 



VI. RESULTS 



We have applied the randomized LP algorithm to test for the different jamming categories 
in practice. The primary aim of this work is not to give exhaustive results, but rather to 
introduce a conceptual framework and some algorithms. Nonetheless, in this section we 
present some sample relevant results for both ordered and disordered periodic packings. 

46 Interior-point algorithms deal better with unboundedness or infeasibility in this context then the simplex 
algorithm. 
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A. Periodic Lattice Packings 



Table 1 in Ref. Jl0[ gives a classification of some common simple lattice packings into 
jamming categories for hard-wall boundary conditions. Table | reproduces this for periodic 
boundary conditions. As we explained in section O, the results in principle will depend on 
the number of unit cells N c chosen as the original packing, and also on the unit cell chosen, 
so the terminology "lattice X is Y jammed" is used loosely here. 



Lattice 


<P 


L 


Z 


C 


N c 


S 




N s 


Honeycomb 


0.605 


Y 


3 


N 


(2,1), (1,2) 


N 


(1,1) 


2 


Kagome 


0.680 


Y 


6 


N 


(1,1) 


N 


(1,1) 


3 


Square 


0.785 


Y 


4 


N 


(2,1) 


N 


(1,1) 


1 


Triangular 


0.907 


Y 


6 


Y 




Y 




1 


Diamond 


0.340 


Y 


4 


N 


(1,1,2) 


N 


(1,1,1) 


2 


SC 


0.524 


Y 


6 


N 


(1,1,2) 


N 


(1,1,1) 


1 


BCC 


0.680 


Y 


8 


N 


(1,1,2) 


N 


(1,1,1) 


1 


FCC 


0.741 


Y 


12 


Y 




Y 




1 


HCP 


0.741 


Y 


12 


Y 




Y 




2 



Table I: Classification of some simple lattices into jamming categories. We give the packing (i.e., 
covering) fraction ip (to three decimal places), the coordination number Z, and the number of 
disks/spheres N s per unit cell, an assessment of whether the lattice is locally (L), collectively (C) 
or strictly (S) jammed (Y is jammed, N is not jammed), and the "smallest" number of unit cells 
N c on which an unjamming motion exists (illustrated at Ref. [21]). 



It turns out that in the cases given in Table |, the packings we have classified as not 
collectively or strictly jammed will not be so for any large N c . Here we give the smallest N c 
for which we have found unjamming motions, and illustrate some of these in Figs. [8] and |9|. 

Moreover, the packings classified as jammed, in this case being the maximal density 
packings in two dimensions (triangular) and three dimensions (FCP and HCP) will be so 
for any finite N c . We leave justification and further discussion of this to Ref. [pT5| . Here 
we just point out for the curious that the triangular lattice is not the only strictly jammed 
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ordered disk packing 47 ; two other examples are shown in Fig. |10| (see Ref. ||). 
B. Periodic Random Packings 

We also tested a sample of periodic random packings in two and three dimensions pro- 



duced via the Lubachevsky-Stillinger compression algorithm [|lT| at different compression 
rates. This algorithm often tends to produce a certain number of rattlers, i.e., spheres 
which are not locally jammed, which we remove before testing for jamming 49 . We would 
like to stress that these are not comprehensive tests, but they do illustrate some essential 
points, and so instead of giving tables with statistics, we give some representative illustra- 



tions. The tolerances for the interparticle gaps (see footnote |22[) used were in the range 
5 G [0.25,0.50], and, as explained earlier, the results in some cases depend on this chosen 
tolerance, but not strongly. 

All random disk (i.e., two-dimensional) packings we tested were not strictly jammed. 
At the typical LS end states of roughly <p « 0.82, we generally found that the packings 
were collectively jammed (with some exceptions such as the packing shown in Fig. |TTJ) , 
although not strictly jammed, as with the packing depicted in Fig. [12|. However, even at 
very high densities (ip ~ 0.89) the packings were only collectively jammed, as illustrated 



in Fig. |13|. Note that quite different properties were observed for the three-dimensional 
packings: All random sphere (i.e., three-dimensional) packings we tested were strictly (and 
thus collectively) jammed. 



47 Conditions for strict jamming and other possibilities for strictly jammed two-dimensional periodic packings 
will be discussed in Ref. ]l8|| , 

48 In the actual LP implementation, we freeze and ignore such particles. 

49 Notice that checking each sphere for local jamming using ([fj]) only once is not enough under this removal 
scheme. Specifically, once a rattling sphere is removed, this may remove some contacts from the packing 
and can make other spheres not locally jammed. Therefore, neighbors of rattlers are recycled on a stack 
of spheres to be checked for local jamming. We have observed that often, particularly in two-dimensional 
LS packings, all disks can eventually be removed on the basis of just the local jamming test starting with 
only a few percent rattlers. 



23 



VII. DISCUSSION 



Our results have important implications for the classification of random disk and sphere 
packings and suggest a number of interesting avenues of inquiry for future investigations. 
Random disk packings are less well-understood then sphere packings. The tendency of 
disk packings to "crystalize" (to form ordered, locally dense domains) at sufficiently high 
densities is well established. For example, Quickenden and Tan |25| experimentally estimated 
the packing fraction of the "random close packed" (RCP) state to be <p « 0.83 and found 
that the packing fraction could be further increased until the maximum value of <p = 0.906 
was achieved for the triangular lattice packing. By contrast, typical random sphere packings 
at ip in the range 0.63 — 0.66 cannot be further densified. 

Our recent understanding of the ill-defined nature of random close packing and of jamming 
categories raises serious questions about previous two-dimensional studies, particularly the 
stability of such packings. Our present study suggests that random disk packings are not 
strictly jammed at (p ~ 0.83; at best they may be collectively jammed. Of course, the old 
concept of the RCP state was invalid in that it did not account for the jamming category 
of the packing. Previous attempts to estimate the packing fraction of the "random loose" 
state p6| are even more problematic, given that this term is even less well-defined then the 
RCP state. The best way to categorize random disk packings is to determine the maximally 
random jammed (MR J) state (see Ref. 0) for each of the three jamming categories (local, 
collective and strict). Such investigations have been initiated [27| and will be carried out in 
the future. 

The identification of the MRJ state for strictly jammed disk packings is an intriguing 
open problem. On the one hand, we have shown that random packings exist with densities in 
the vicinity of the maximum possible value (ip = 7t/2a/3) that are not strictly jammed, and 
on the other hand, there is a conjectured achievable lower bound cp > ^/3n/8 corresponding 
to the "reinforced" Kagome lattice (see Fig. pi]). For random sphere packings, an initial 
study undertaken in Ref. [p7|| , using the LP algorithm described in this work, found that 
maximally disordered random packings around ip ~ 0.63 were strictly jammed. This suggests 
a close relation between the conventionally accepted RCP packing fraction and the packing 
fraction of the MRJ state for strictly jammed packings. However, it has been shown |27[ that 
at a fixed packing fraction ip « 0.63 the variation in the order can be substantial and hence 
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packing fraction alone cannot completely characterize a random packing. The conventionally 
accepted RCP packing fraction in two dimensions may be approximately close in value to 
the packing fraction of the MRJ state for collectively jammed packings. Much less obvious 
is what is the MRJ state for collectively jammed sphere packings. Finally, a completely 
unexplored question concerns the identification of the MRJ state for locally jammed disk 
and sphere packings. 

VIII. CONCLUSIONS 

In this work we have proposed, implemented, and tested a practical algorithm for verifying 
jamming categories in finite sphere packings based on linear programming, demonstrated 
its simplicity and utility, and presented some representative results for ordered lattices and 
random packings. Interestingly, the random packings that we tested were strictly jammed 
in three dimensions, but not in two dimensions. Future applications of the randomized 
linear-programming algorithm are to be expected. We will further present and explore the 
theoretical connections between rigidity and jamming, kinematic and static rigidity, rigidity 
and energy, rigidity and stability, and finite, periodic and infinite packings in Ref . |H| , and 
work is already under way to provide highly efficient implementations of various optimization 
algorithms for linear and nonlinear programming on large-scale (contact) networks. 
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IX. FIGURE CAPTIONS 



Figure 1: Unjamming the honeycomb lattice inside a hard- wall simple box container. The arrows 
in the figures given here show the direction of motion of the spheres V in the linear unjamming 
motion, scaled by some arbitrary constant to enhance the figure. N c = (3, 2) unit cells make this 
small packing. 

Figure 2: Unjamming the honeycomb lattice. A sub-packing of size N c = (3, 3) of an infinite 
honeycomb lattice packing is pinned by freezing all neighboring image disks. A representative 
unjamming motion is shown as a sequence of several frames between times t = and t = 1. The 
unshaded disks represent the particles in the generating packing -P(R), while the shaded ones are 
image disks that touch one of the original disks. 

Figure 3: Unjamming the Kagome lattice. Periodic boundary conditions are used with N c = (2, 2). 

Figure 4: Feasible displacements polyhedron. The figures show three stationary (dark gray) disks 
surrounding a mobile disk (light gray). For each of the three stationary disks, we have a nonlinear 
impenetrability constraint that excludes the mobile disk from a disk of radius D surrounding each 
stationary disk (dark circles). Also shown are the linearized versions of these constraints (dark 
lines), which are simply tangents to the circles at the point of closest approach, as well as the 
region of feasible displacements bounded by these lines (shaded gray). 

This region is a polyhedral set, and in the left figure it is bounded, meaning that within the ASD 
the mobile disk is locally jammed (trapped) by its three neighbors, while on the left it is 
unbounded, showing the cone of locally unjamming motions (escape routes). Notice that with the 
true nonlinear constraints, the mobile disk can escape the cage of neighbors in both cases, 
showing that the ASD is not exact. However, it should also be clear that this is because we have 
relatively large interparticle gaps here. 
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Figure 5: The packing from Fig. |4| shown again with a numbering of the disks. D, L denotes particle 
i and denotes the contact between the ith and jth particles, i.e., the contact {i, j}. 

Figure 6: Contact network of a random packing of 100 disks with periodic boundary conditions 
and 5 = 0.5. Periodic contacts with neighboring image spheres are also shown. All disks are locally 
jammed within the rather large gap tolerance employed. 

Figure 7: Example of a lattice deformation. The above periodic packing (packing 3 in Ref. |22[| ) is 
collectively jammed, but not strictly jammed. It can be continuously sheared toward the triangular 
lattice by deforming the lattice in a volume-preserving manner, as shown here. 

Figure 8: Simple collective mechanisms in the Kagome and honeycomb lattices, respectively. These 
lattices are not collectively jammed with periodic boundary conditions, as the sample unjamming 
motions periodic with N c = (1, 1) for Kagome (left) and N c = (1, 2) for honeycomb (right) shown 
here illustrate. 

Figure 9: Shearing the honeycomb lattice. The honeycomb lattice is not strictly (or collectively) 
jammed, and an example of a lattice deformation with N c = (1, 1) is shown, replicated on several 
unit cells to illustrate the shear character of the strain e = (AA)A- 1 [c.f. (0)]. Note that only 
three (original) spheres are involved in the actual calculation of this unjamming motion, the rest 
are image spheres. 

Figure 10: Examples of strictly jammed lattices in two dimensions. The 6/7th lattice (packing 
number 2 in Ref. [^] and the last packing in Ref. [||), left, is obtained by removing every 7th 
disk from the triangular lattice. The reinforced Kagome lattice, right, is obtained by adding an 
extra "row" and "column" of disks to the Kagome lattice and thus has the same density in the 
thermodynamic limit, namely, it has every 4th disk removed from the triangular packing (see also 
Ref. @). 

Figure 11: A random packing {ip = 0.82) of 1000 disks that is not collectively jammed, and a 
representative periodic unjamming motion. 
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Figure 12: A random packing (</? = 0.83) of 1000 disks that is collectively jammed but not strictly 
jammed, and a representative unjamming motion. Though it is hard to see from this figure, this 
is indeed a shearing motion that induces an unjamming mechanisms. A more insightful animation 
can be found at the webpage [21]. 

Figure 13: A dense (ip = 0.89) random packing of 1000 disks that is collectively jammed but not 
strictly jammed, and a representative unjamming motion. One can see the grains gliding over each 
grain boundary due to the shear, bringing this packing closer to a triangular lattice. 

X. FIGURES 
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Figure 1: Donev et al. 
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Figure 2: Donev et al. 
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Figure 3: Donev et al. 
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Figure 5: Donev et al. 
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Figure 6: Donev et al. 
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Figure 7: Donev et al. 
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Figure 8: Donev et 
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Figure 10: Donev et al. 
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Figure 11: Donev et al. 
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Figure 12: Donev et al. 



41 




Figure 13: Donev et al. 
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